Examples of Synchronization in Discrete Chaotic Systems 



Juan C. Botero and Jean- Jacques E. Slotine 



Abstract — This paper presents an application of partial 
contraction analysis to the study of global synchronization 
in discrete chaotic systems. Explicit sufficient conditions on 
the coupling strength of networks of discrete oscillators are 
derived. Numerical examples and applications to simple systems 
are presented. Previous researches have shown numerically 
that the systems under study, when arranged in a network, 
exhibits rich and complex patterns that can dynamically change 
in response to variations in the environment. We show how 
this "adaptation" process strongly depends on the coupling 
characteristics of the network. Other potential applications of 
synchronized chaotic oscillators are discussed. 

I. INTRODUCTION 

This paper presents a research focused on coupled discrete 
dynamical systems [1], [9]. In particular, it is of interest 
to study the synchronization phenomenon in chaotic sys- 
tems, i.e. systems in which small differences on the initial 
conditions can lead to a completely different behaviours 
in time. In spite of this characteristic, synchronization can 
be achieved in this type of systems [10], [11]. This work 
presents a number of examples in which synchronization 
can be guaranteed by choosing the appropriate coupling 
strength. The coupling strength can be determinate using 
partial contraction analysis. This analysis is able to pre- 
dict complete synchronization independently of the initial 
conditions, as long as the coupled system verifies certain 
properties. Besides synchronization, diffusively coupled dy- 
namical systems oscillating chaotically in time have seen to 
lead to interesting emerging properties due to changes in the 
environment [5], This feature has been used in this work to 
simulate the motion of a two-legged robot where each leg 
is a chaotic oscillator and these oscillators are diffusively 
coupled. The appeal of such a system is the fast adaptation 
it shows when sudden unexpected dynamic changes in the 
environment occur. This method of locomotion does not 
require any control or optimization process to be performed. 
Furthermore, this type of analysis can be done on any 
periodical task that requires fast responses to unpredictable 
changes. 

The role of chaos in "adaptation" and evolution has been 
widely discussed before [4], [8]. In particular, some research 
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in this area, sometimes known as "adaptation at the edge of 
chaos", suggests the possibility that when biological systems 
adapt in order to survive, the process of evolution may favor 
those systems that are near a phase transition from order 
to chaos [8]. In any event, the chaotic behaviour of simple 
dynamical systems can be exploited to obtain adaptation 
behaviors to external dynamical disturbances. 

The method used is mathematically simple and also similar 
to for the analysis of coupled dynamical systems [3], [9], 
[10]. Defining proper matrices to project the states onto 
the so-called synchronization manifold [10], [11], contrac- 
tion theory [7], [13], [14], [17] can be used to determine 
conditions for which a set of coupled identical systems will 
completely synchronize. 

II. Mathematical Foundations 

As in [13], consider a set of coupled dynamical systems 
and define the vector x{}: 
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so < Uj.Vj > = 0, where <x,y> denotes standard inner 



product. Thus, the elements of the vectors x" 
orthogonal "one to one": 
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III. Contraction Theory for Discrete Systems 

Some basic results in Contraction Theory are presented as 
follows: 



Definition 1 

Given the systems of equations x i+1 = fa (xt,i) , a region 
of the state space is called contraction region with respect 
to a uniformly positive definite metric Mi (x i} i) = if 
in that region [7]: 
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where 



f = e t 
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Theorem 1 



Given the systems of equations x i+i = fi(xi,i) , any 
trajectory, which starts in a ball of constant radius with 
respect to the metric Mi(xi,i), centered about a given 
trajectory and contained at all times in a contraction region 
with respect to the metric Mi(xi,i), remains in that ball 
and converges exponentially to this trajectory. Furthermore, 
global exponential convergence to the given trajectory is 
guaranteed if the whole state space is a contraction region 
with respect to the metric Mi(xi 1 i). This corresponds to 
a necessary and sufficient condition for exponential conver- 
gence of the system [7]. 

Theorem 2 

Consider two coupled systems. If the dynamics equations 
verify 

- h (x^) = - h (*<?>) (8) 

where the function h is contracting, then x^ and x^ 
will converge to each other exponentially, regardless of the 
initial conditions [17]. 

Hence, the process of one variable converging exponen- 
tially to another, i.e. synchronization, can be seen equiv- 
alently as contraction of the projection in the transverse 
manifold [13]. It is clear since the transverse manifold is, by 
construction, perpendicular to the synchronization manifold: 



xi — > x 2 <S=> x 







(9) 



If the projection of a particular solution has null norm in 
the transverse manifold, it lays entirely on the synchroniza- 
tion manifold (here there are no distinction between synchro- 
nization and oscillator dead). Thus, it is only necessary to 



study the hereafter called error dynamics Xi — Xj, and prove 
its contracting behaviour. The proof can be seen in [13] and 
it can be express in rather intuitive terms: synchronization 
will be achieved if the difference of the states of any two 
oscillators is zero. 



IV. Application of Contraction to Chaotic 
Systems Synchronization 

A. General Case: Chaotic Maps 

In order to study the phenomenon of complete synchro- 
nization in chaotic systems, it is usual to work with one- 
dimensional maps [5], [12] due to the complex behaviours 
that can be observed from such simple models. When no 
coupling is present between any of the maps, due to their 
chaotic nature the systems behave random-like. But when 
appropriate coupling is build between the maps, one will 
expect some sort of interaction. Hence, in the case of 
synchronization or complete synchronization, it is of interest 
to have a coupling with two main properties [12]: 

1) The coupling must make the states of the systems 
closer to each other, i.e. dissipative coupling. 

2) The coupling must not affect the synchronization state. 
A general form of this kind of coupling operator of 

dimension-2 is: 



1 - a (3 
(3 1-/3 



(10) 



where < a < 1 and < b < 1. The simplest case 
of this coupling would be a = (3 = e [12]. This case 
allows the coupling to be symmetric, which leads to a great 
simplification in the computations needed, and furthermore, 
it allows the direct application of the proposed method. 
Hence, it is possible to couple any two maps x and y as 
[12]: 
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and with a proper choice of the coupling strength e, 
the complete synchronization state can be achieved. For 
example, two oscillators following the well known skew map 
[8], [12]: 

f 2il if o < x n < a 
/<*.) = { *"„<-<! 

which exhibits pure chaotic behaviour due to the subse- 
quent stretching and folding processes in the interval from 
to 1 [12]. Two discrete systems of this type, say x and 
y, can be coupled as in (11). Thus, varying the coupling 
strength e it is possible to achieve complete synchronization 
of the states. The complete proof of this example can be 
found in [12], where the well known technique of transverse 
Lyapunov exponents is used [10], [11]. 



However, theorem 2 from contraction theory provides a 
simple yet powerful tool to analyze this particular system. 
The following theorem was presented for the first time in 
[17]. Using this theorem it is possible to write the system of 
two skew maps x and y, coupled together in the following 
form: 

x n+1 - (1 - 2e) / (x n ) = y n+l - (1 - 2e) f (y n ) (12) 
Hence, it is only necessary to verify that the function 




Fig. 1. Chaotic oscillators coupled with different strengths. Left: e — 0.9, 
i.e. outside of the contraction range. Right: e = 0.4, i.e inside the predicted 
synchronization interval. 



h (x) = (1 - 2e) / (X) 



(13) 



is contracting for f(x) in (12) (for a particular metric). 
As df(x)/dx is not continuous, two cases must be taken into 
account. Each one of these cases should give an interval for 
e in which contraction of the function h(x) is guaranteed. 
The general result is the intersection of these two intervals 
(in fact, the supremum or the biggest set contained in the 
intersection of the intervals obtained). In this particular 
example, we must show that the Jacobian 



F = 



dh 0) 
dx 



(1-20 
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dx 



B. Coupled Map Lattice ( CML) 

In [5], the authors proposed a more general model a 
coupling in which the states takes into account the dynamics 
of any other oscillator "next to it" equally weighted: 



4+1 = (!-£)/ (4) + 2 (/ (4 +1 ) + / (4^)) d5) 

As it can be seen, at the iteration n+1 the dynamics of the 
i-th oscillator depends on its previous state (weighted by 1-e) 
and the previous state of the "nearest" two oscillators, e.g. (i- 
Ij-th oscillator and the (i+l)-th oscillator, equally weighted 
(14) by a factor of e /2. Let x n +i= f(x) =1 - ax\ where a is a 
positive constant and write the system in the following form: 



is contracting in some metric M. Indeed, if M is equal 
the identity matrix /, for the two cases we find the suitable 
intervals 



< x < a 
a < x < 1 



df(x) 
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Thus, the supremum of the intersection of the two intervals 
is the second one (as it is contained in the first interval). It 
is straight forward to verify that the use of the transverse 
manifold leads to the exactly same conditions as before, i.e. 
an interval for e in which the error dynamics is contracting 
in the same metric as before. In fact, the two synchronization 
intervals obtained for the error are 



< x < a 



a < x < 1 
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< £ < 1 



Going beyond the analysis done by transverse Lyapunov 
exponents, the contraction theory is capable of giving in- 
tervals where the selected coupling strength e guarantees 
complete synchronization. This can be clearly seen in figure 
1. For two coupled oscillators X and Y and a = 0.7, 
the proposed method predicts complete synchronization (in 
which the values lie on the X = Y line) for a coupling strength 
between 0.35 and 0.65. Fig. 1 shows two cases. 



4+i = / (4) + 2 (/ (4 +1 ) + / (4 _1 ) - 2/ (4)) 

(16) 

Thus, for two coupled oscillators x and y, the matrices U 
and V are vectors 

U = [ 1 1 ] 
V = [ 1 -1 ] 

and the dynamics of the projection on the synchronization 
manifold x = e = x - y is given by 



<Sn+l — Xn+1 — lln+1 
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(17) 



where a = x+y, and no approximations, e.g. linearization 
about a certain point, have been made. Now, a virtual 
auxiliary system can be constructed [17] 
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where both and e„ are particular solutions. This new map 
shall be contracting, i.e. df{w n )ldw n is u.n.d , in any case 
for an identity metric, i.e. in (7) the matrix M = I (identity 
matrix), and therefore e will tend exponentially to if [7]: 



2(a-A) 2 (a+i) 

< £ < 
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Letting a to be variable, (19) describes a family of 
hyperbolas in the a - e plane. Now, it is necessary to find a 
region in the a - e plane that guarantees contraction, i.e. a 
region where any pair (a - e) will assure contraction of the 



system. This can be easily done by finding the supremum of It is worth noticing that this is the same result on would 



the areas defined between the curves of the hyperbolas. In 
Fig. 2 some of the hyperbolas are shown for a = 1.7, where 
the map exhibits a chaotic behaviour (a± < a% < 03). 




Fig. 2. Family of hyperbolas defining contraction regions for the coupled 
system. The hyperbola with the "narrowest" area between the curves (0-3) 
defines the sufficient contraction region for the error dynamics 



find applying theorem 2: writing the system of equations in 
the form: 
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Vn+l 
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(22) 



where, in order to assure synchronization of the states, the 
function: 



(y 



1 



(23) 



must be contracting in some metric. Using the identity 
metric the same results as before are found. 

Thus, the expected behaviour should be (looking for 
positive values of e ) as presented in figure 4, where the 
shaded area represents the guaranteed contraction region 
for the coupled system (recall that each oscillator without 
coupling follows its own nominal chaotic dynamics): 



In this particular case, it is found that the supremum is 
given by the area between the curves o the hyperbola with 
a =2(2/^/a). Now, since the initial value of x is random 
number between and 1, the absolute values of x and y 
are always less than (y/2/a) « 1 (for values of a which 
guarantee chaotic behaviour), it can be easily shown that: 



-1 < 1 — ax < 1 => -2 < —ax < 



(20) 



Since both a and x 2 are always positive, the condition in 
(20) is true. 

The values of a which guarantee chaotic behaviour can 
be found by computing the Lyapunov exponents of the map 
as a varies [16]. The first Lyapunov exponent for different 
values of a is shown in Fig. 3. 
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Fig. 3. First Lyapunov exponent of the logistic map as a varies 

It can be seen that the values close to 2 guarantee chaotic 
behaviour. Hence, the limit, i.e. the narrowest hyperbola is 
given by a w 2 , since 



2a - 1 2a + 1 

< e < 



3 a 



3 a 



(21) 




Fig. 4. The shaded area is the guaranteed contraction region for a- 2 and 
an identity metric 

Numerical experiments were performed: using a Sobolev 
distribution (which gives a random but uniformly distributed 
cloud of points in a determined region) in the a-e plane, Fig. 
4 shows the norm of the error vector for each couple (a,e) 
after a transient time of 300 iterations. We find indeed, as 
shown in figure 5, that the actual region of guaranteed con- 
traction for the coupled system resembles the one predicted 
by the proposed method. 




Fig. 5. Actual contraction region for a = 2. The dark area represents a set 
of pairs (a,e ) for which the norm of the error is zero after the transient 



V. Simulations on the Coupled Map Lattice 

Using the previous results, the coupling strength can be 
tuned to guaranteed synchronization, as shown in Fig. 6. 
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Fig. 6. Two chaotic discrete oscillators. Above: non coupled. Below: 
properly coupled 



From these results, it is clear how, using a value for 
the coupling strength that satisfy the conditions for the 
contracting error dynamics in the synchronization manifold, 
the two systems will tend exponentially to each other (equiv- 
alently, the error will tend exponentially to the origin of the 
synchronization manifold). 
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Fig. 7. Outline of the model. This model is inspired on the one proposed 
in [5] 




Fig. 8. Model used for the two-legged robot simulation 



VI. Adaptation and Emerging Behaviour Via 
Sensing 

Let's modify the coupled map lattice considered above 
such as the system receives information of the environment 
via sensing. The signal coming from the sensors can be any 
kind of signal (rotation angle, voltage, etc) as long as it can 
be related somehow to the state variable x. Let's consider for 
example the following dynamics as proposed in [5]: 



+ E1 



(24) 



where s z n denotes the sensor reading at the iteration n and 
s„ is the mean value of the sensors at the iteration n. The 
meaning of this type of coupling is as follows: each oscillator 
follows its own the chaotic dynamics, but this is "adjusted 
and updated" to reduce the global (via £i) and local (via 
£2) difference in the sensors [5]. Fig. 7 shows the model 
schematically. 

Let's suppose that each oscillator represents the dynamics 
of a motor and that each motor is driving a leg of a 
simple two-legged robot. Thus, each leg has an actuator 
following the proposed dynamics coupled in order to achieve 
synchronization. Let the variable x be the command to the 
actuator (torque r) and the signal sensed (s) be the angle 6. 
This to variables are related as follows in Fig. 8: 



Following the method as done previously in (18), the error 
dynamics is given by: 
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For this is particular case, the angle s can be related 
directly to the torque x in a linear fashion, e.g. in (25) 
s = mx + b (where the constants m and b are given by the 
sensor characteristics). Hence once again, a virtual auxiliary 
system can be construct 

Vn+i = f (y n ) = -aay n + 2a (y + e 2 ) a(s^ - s^) 

where and e n are particular solutions. Using an identity 
metric M = I, the system shall be contracting, i.e. the error 
e will converge exponentially to 0, for: 
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(26) 



i> = (K* -Nl + mgL) 



and K* is the equivalent stiffness of the two springs. 
Fig. 9 shows two coupled oscillators starting at random 
initial conditions. The complete synchronization is achieved 
by tuning the coupling strength to be inside the predicted 
synchronization range. 



n n 

Fig. 9. Proper selection of the coupling strength assures complete syn- 
chronization. Left: coupling strength outside the contraction range. Right: 
coupling strength inside the contraction range 



x[m] 

Fig. 11. Friction coefficient between leg and surface decrease after y = 
1. Left: decrease of 25% in the friction coefficient in the left leg. Right: 
decrease of 40% in the friction coefficient in the right leg 



VII. Why Bother with Chaotic Systems? 

At this point it seems natural to ask why it is worth 
working in a chaotic regime, since similar results, specif- 
ically synchronization, can be achieved with non chaotic 
oscillators. Apparently the key to answer this questions is 
adaptability. The following example shows the advantages of 
a chaotic regime. Assume the two-legged robot under three 
different conditions: first assume that the legs are chaotic but 
not coupled; second, the legs are coupled but not chaotic; and 
finally the legs are coupled working in the chaotic regime. 
As an external condition assume that the left leg senses a 
decrease of 25% in the friction coefficient after y = 1 m. 
Figure 10 shows the results obtained by simulation. 



- 



y [m] 

1.5 







O XsCS;a-jzic m Coupling 
x Chaotic m ; Coupling 
□ Chaotic wa Coupling 


























i i 


]=. i i i 



-1.5 -J -0 5 5 1 



x [m] 

Fig. 10. Trajectory followed by the robot in the x-y plane 

Apparently, the condition coupled-chaotic is more sensi- 
tive to external changes in the environment, which can be 
seen as an emerging adaptive property of the system. 

Another simulation under the condition coupled-chaotic 
was done in which after some fixed traveling distance one 
of the legs "feels" a change in the surface, e.g. in this case 
the sensor detects a change in the angle due to a decrease in 
the friction coefficient between the leg and the surface. The 
results suggest that the change in the behaviour depends on 
the magnitude of the change in the environment. 

These results suggest that the idea of adaptation to changes 



in the environment is intrinsic to the system dynamics. In the 
sense of adaptation and evolution [5] and perhaps abusing 
of the terminology, chaotic oscillations may be regarded as 
"mutants" from a nominal (non chaotic) oscillation pattern. 
Such deviations from the smooth nominal behaviour give 
the system the capacity of response faster to external inputs. 
Furthermore, due to the nature of chaotic oscillations, e.g. 
strange attractors where the trajectories are random-like yet 
confined in space [15], they appear to be a good compromise 
between a "mutation" sufficiently strong to make the system 
adapt to a new external condition but mild enough to 
maintain information from the previous behaviour. 

Finally, note that the analysis does not require any goal 
direction, because the system, due to its simple configuration, 
tends to "walk straight ahead". This feature makes unnec- 
essary the implementation of any mean of command, e.g. 
algorithms deciding to freeze a leg at a certain moment. The 
system simply moves forward, responding quickly to changes 
in the environment and, more interesting, this response is 
proportional to the magnitude of the change. 

But the reason why the chaotic regime exhibits this 
adaptive behaviour remains an open question. In order to 
understand this phenomenon let's take the following exam- 
ple: assume that each movement of the legs takes place in 

1 second (arbitrary time interval where one oscillation takes 
place). It is obvious that the power spectral density of the 
signal (PSD) will have, in any case, a peak at 0.5 Hz (it takes 

2 seconds for the leg to complete a cycle forward-backward). 
Fig. 12 shows the mean PSD of 128 runs starting at 

random initial conditions: the mean of 128 non chaotic 
signals, i.e. with a value of a that makes the first Lyapunov 
exponent negative, e.g. a = 1.1 (dashed line) and the mean 
of 128 chaotic signals, i.e. a = 1.7 (solid line). Both signals 
are subjected to a equivalent coupling strength 

It can be seen how the chaotic signal carries energy at 
some other different frequencies than the trivial one. In par- 
ticular, it can be seen that the period-doubling phenomenon 
occurs (peak at half the main frequency). Due to its chaotic 
nature, this particular feature allows the system to have 
richer movement patterns in changing environments (because 
different frequencies are being excited by the same signal). 

For example, consider the case of a small vehicle where 
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Fig. 12. Power spectral density of a chaotic and a non-chaotic oscillator 



each wheel is driven by a chaotic oscillator and these oscil- 
lators are coupled as before. As any road profile irregularity 
can be describe as an ergodic random process in space [2], a 
vehicle traveling at constant speed V (in average) will sense 
an input which is an ergodic random process in time. As the 
energy of the signal in chaotic regime is amplified at more 
particular frequencies, the natural response of the system, 
e.g. the legged robot or the wheeled vehicle, will show richer 
outputs. This could be seen as an advantage when dealing 
with challenging and changing environments. 

VIII. Networks of Symmetric Coupled 
Oscillators 

Up to this point, this work has focused on two coupled 
oscillators. However, the extension to more oscillators in a 
symmetric network is straight forward. Indeed, as it will 
be shown in this section, the synchronization analysis of 
this particular type of networks can be done as if it were 
composed by two oscillators. 

Let's start by saying that matrices U and V in (2) and (3) 
can always be defined for n identical oscillators with scalar 
dynamics and symmetric coupling, leading to a synchroniza- 
tion of the states of any two adjacent oscillators. 

Now, suppose n identical oscillators, coupled in such 
way that the each oscillator updates its own dynamics by 
sensing the dynamics of every other oscillator in the network. 
Furthermore, suppose that the magnitude of the coupling 
strength is equal for any pair of oscillators in the network. 
The structure of a network with verifying such conditions is 
known as "all-to-all" symmetry structure [17]. As a general 
case, consider the Coupled Map Lattice: 
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where f(x\ ) is any function representing a discrete 
dynamical system (for k = 1, 2,. . . , n) and g(e) is a general 



function of the coupling strength. Using (0), it is possible to 
find the projection of the system on the transverse manifold, 
leading to a (n — 1) dynamics of the form 



-l = Te„ 



(28) 



where T is a diagonal matrix where the diagonal elements 
are functions of the form Vjj(s, x% , xjf' ) for j = 1,2,..., 
n-1. Thus, each element of the diagonal is a function that 
can be expressed for any pair of oscillators in a "two-by- 
two" fashion. If Vjj is upper bounded, then the coupling 
strength e can always be tuned to guaranteed uniformly 
negative definiteness of the diagonal matrix T [17]. It is 
worth noticing that in the case of chaotic oscillators the 
condition on the upper bounded of the function is guaranteed 
by the strange attractor behaviour [16], as agreed in section 
VIII. 

This results leads to a convenient simplification in the 
analysis of any symmetric network of oscillators, since it is 
necessary to study one single function g to guarantee com- 
plete synchronization. Furthermore, it can be stated that the 
synchronization analysis of any network of oscillators with 
an "all-to-all" symmetry structure can be made studying the 
uniformly negative definiteness of a single scalar function. 
For example, for three discrete oscillators x, y and z: 



in+i = / (a;) = 1 
Vn + l = f (v) = 1 ' 
Zn+l = / (z) = 1 - 
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working in chaotic regime, e.g. a ■ 
operator in the form 
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can be analyzed as previously, e.g. x± = = x - y and 



= (2) 



r e (l) " 




e (2 > 

. n + l . 





y — z leading to 



(e — 1) acri 

(e-1) ao-2 



-"XX 

„(2) 



(31) 



where a\ = (x + y) and 02 = (y + z). Once again, a virtual 
auxiliary system in the form 



[ Zll ] - [ 
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for which the vectors and e„ are particular solutions. 
Thus, the states will converge exponentially to each other 
for a coupling strength e that guarantees uniformly negative 
definiteness of the matrix F, i.e. the largest singular value of 
F remains smaller than 1 uniformly [3]. This is easily done 
(due to the diagonal form of the matrix F) by tuning the 
parameter e in the function 



7(e) = (e- l)a k a 



(33) 



so that 7 is uniformly negative definite. The upper bounded 
value of (Tfe is easily checked from the previous examples 

(<7fc~2). 

Finally, for more general systems, the system of equations 
in d29t can always be expressed in such way that the 
following theorem from classical contraction analysis can be 
directly applied [17]. 

Theorem 3 

Consider q coupled systems. IF a contracting function 
h{xi) exists such that 

=^l- h {^) (*<?>) 

THEN all the systems synchronize exponentially, regard- 
less of the initial conditions. 

In the case of the last example, the system dynamics can 
be expressed as 

x t+ i - (e - 1) ax\ = yt+i - (e - 1) ayl = z t+1 - (e - 1) az^ (34) 

for which the analysis of the function h(xt) = (e - l)aa;| 
and the tuning of e to achieve uniformly negative definiteness 
of h, leads to the same result as above. 

IX. CONCLUSIONS AND FUTURE WORKS 

A. Conclusions 

This paper has shown the feasibility of studying synchro- 
nized discrete chaotic oscillators, by means of a differential 
analysis of the nonlinear system. Different kinds of networks 
and coupling operators can be analyzed globally by simple 
method. Further investigations on the emerging adaptation 
properties of such systems are necessary for practical imple- 
mentations. 

B. Future Works 

Future work will be focus on the extension to systems 
with more degrees of freedom. Some of the work done on 
synchronization is intended to be implemented in a set of 
Lego Mindstorm® . A small two-wheel vehicle was built 
using the set of Lego® blocks and the input signal to the two 
motors driving the wheels is given by (26). The algorithm 
is intended to be implemented using the angle of rotation as 
the sensed signal, i.e in accordance with the model of the 
two-legged robot presented in simulations. 
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